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Abstract 

The observation of a small change in spectral slope, or 'knee' in the fluxes of cosmic 
rays near energies lO^^eV has caused much speculation since its discovery over 
40 years ago. The origin of this feature remains unknown. A small workshop to 
review some modern experimental measurements of this region was held at the 
Adler Planetarium in Chicago, USA in June 2000. This paper summarizes the results 
presented at this workshop and the discussion of their interpretation in the context 
of hadronic models of atmospheric airshowers. 
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1 Introduction 

The primary cosmic ray particles extend over at least twelve decades of energy 
with a corresponding decline in intensity of over thirty decades. The spectrum 
is remarkably featureless with little deviation from a constant power law across 
this large energy range. The small change in slope, from oc to oc E~^'^, 

near 10^^ eV is known as the 'knee' of the spectrum. The coincidence of this 
feature with the highest energies expected from diffusive shock acceleration in 
supernova remnants is intriguing [1]. It can also be attributed to contributions 
from nearby or recent supernova events [2]. To distinguish between these var- 
ious ideas, better measurements of the elemental composition of cosmic rays 
at these energies are essential. The accurate determination of composition has 
provided some of the key advances in this field at lower energies, where direct 
measurements are possible with detectors above the atmosphere. For example, 
the realization that the observed spectral slope of cosmic rays is significantly 
steeper than that produced in the cosmic ray sources themselves resulted from 
measurements with sufficient elemental resolution to separate primary source 
cosmic ray elements from those produced in the interstellar medium at 10^^ 



Here we collect together various recent measurements of cosmic rays in the 
'knee' region which were discussed at a workshop held in Chicago in June 



an overview of the state of experimental measurements in this field and to 
provide some suggestions for future improvements. The systematic limitations 
of the type of measurements presented here are also evaluated. In the modern 
era, these limitations are associated in part with the validity of the various 
numerical simulations used to interpret the data and the absolute energy scale 
for the data. 

The structure of this paper basically keeps to the structure of the workshop 
where both experimental results and various aspects of the numerical simu- 
lations were discussed. This work is not an exhaustive summary of the field 
but summarizes the efforts that were represented at the workshop. There are 
other notable experiments in this area which were not discussed in the work- 
shop, in particular HEGRA[4] and EAS-T0P[5]. Also, developments since the 
workshop are not discussed here, these are summarized in the proceedings of 
the 27th International Cosmic Ray Conference (2001) Hamburg[6]. 

We begin with a brief overview of the relative merits of the experimental 



eV [3]. 




paper is to provide 
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Fig. 1. Participants in the cosmic ray workshop at Adler Planetarium. (Left to 
Right) Moissev, Matthews, Swordy, Fortson, Sills, Cherry, Minnick, Horandel, Kam- 
pert, Coutu, Knapp, Roth, Kieda, Wakely, Pryke, Miiller, Hinton, Fowler, Cronin, 
Spiczak, Kettler, Shibata, Cao 

techniques followed by short descriptions of direct and indirect (air shower) 
measurements presented at the meeting. There follows a relatively detailed 
discussion of air shower simulations and the physics which goes into them. We 
then present a simple comparison of the results of the experiments in terms of 
the mean logarithm of the inferred primary mass. Contemporary experiments 
are breaking new ground in terms of statistics and modeling, and we include 
some highlights of multi-species fits and multi-parameter correlation studies. 
Finally we draw some conclusions. 



2 Overview of Techniques 



This section seeks to summarize the main cosmic ray detection techniques 
and their relative merits. By no means comprehensive, the points that are 
touched on here were felt to be most relevant to the experiments discussed at 
the meeting. The experiment acronyms are defined in Section 3. 
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2.1 Direct Detection 



The detection of cosmic rays above the atmosphere is the only way to obtain 
direct measurements of the primary particles and their energy spectra. There 
are a number of techniques for direct detection. So-called "active" techniques 
usually involve some combination of charge detection instruments and energy 
measurement instruments, such as scintillator combined with a calorimeter or 
transition radiation detector (TRD). Here, the data can be recorded electron- 
ically in-flight. TRD's rely on the passage of a particle through the detector 
without a nuclear interaction while calorimeters not only require an interaction 
but must contain as much of the ensuing cascade as possible for accurate en- 
ergy measurements. It should be noted that present TRD's saturate at about 
10^^ eV and provide better energy measurements of high Z nuclei compared to 
protons and Helium. In contrast, calorimetric techniques can be used to obtain 
reasonable energy measurements of protons and Helium up to knee energies. 

Nuclear emulsion is an example of a passive technique that, when interleaved 
with suitable converters, can simultaneously measure the primary charge and 
energy. These payloads must be recovered and subjected to complicated and 
exacting procedures to extract the energy and charge information. While each 
technique has its own strengths and weaknesses, the major limitation to all 
direct techniques to date is the inability to fly large area detectors for long pe- 
riods of time. This limit is constantly being pushed as with the recent advent 
of Ultra-Long-Duration-Balloon flights and by constructing payloads that do 
not require pressurized vessels for operation. However, the nuclear emulsion 
experiments currently have the largest exposure and therefore the most sta- 
tistically significant data at the highest energies. It is these experiments that 
will be discussed below. 



2.2 Indirect Detection 



Once the primary particle interacts with the Earth's atmosphere, the detection 
techniques indirectly measure the primary parameters through the air shower 
components produced in the interaction. It is then a matter of reconstructing 
the primary energy and nuclear charge from the air shower observables. The 
main difficulty with indirect techniques is that these reconstructions depend 
on hadronic interactions which are empirically undetermined at this time for 
the relevant energies and kinematical regions. However, large area, long du- 
ration experiments can be operated from the ground allowing statistically 
significant measurements to be made up to the highest energies in the cosmic 
ray spectrum. 
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The different techniques can be separated into detectors that measure the par- 
ticle content of the shower (KASCADE, CASA-MIA), and those that measure 
Cherenkov hght (DICE, BLANCA, VULCAN). 

The slope and normalization of the lateral distribution of an air shower com- 
ponent are functions of the energy and mass of the primary cosmic ray. While 
composition information can be obtained by lateral distribution measurements 
from a charged particle array (e.g. the CASA-MIA results below use the elec- 
tron lateral distribution slope), in practice, it is much easier to make density 
measurements of the required precision when working with Cherenkov light. 
Because many photons are emitted by a given charged particle, and these 
photons are not significantly attenuated by the atmosphere, they are much 
more numerous than charged particles at ground level. Additionally, since the 
majority of the light is emitted well above the ground, and at significant an- 
gles to the shower axis, its lateral distribution at the ground is flatter than 
that of charged particles. For these reasons, the sensitivity and dynamic range 
requirements of Cherenkov detectors are not as stringent as those for charged 
particle detectors to probe the primary composition in a given energy range. 
These are the strengths behind the BLANCA and VULCAN results discussed 
below. However, this advantage is at least partially offset by the requirement 
that Cherenkov detectors can only run on dark, moonless nights. 

The DICE experiment measures the angular distribution of Cherenkov light, 
and thus makes a somewhat direct measurement of the longitudinal profile of 
the air shower. This puts it a different class to all of the other experiments 
discussed hertPI which measure quantities which have been "integrated" over 
the shower development history. 

The multi-component nature of air showers allows event-by-event cross-correlation 
measurements. This approach has been widely exploited by the KASCADE 
experiment, it is also possible with the DICE and BLANCA experiments which 
are co-located with CASA-MIA (see Section 5.3). 



3 Overview of Experiments 

Here we provide a brief description of the measurements discussed at the 
meeting, including direct detections by RUNJOB and JACEE and the air 
shower experiments KASCADE, CASA-MIA, DICE, BLANCA and SPASE- 
VULCAN. 



With the exception of the KASCADE low energy muon tracking detectors which 
are not mentioned in this paper. 
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3.1 JACEE 



The Japanese-American Cooperative Emulsion Experiment (JACEE) is a se- 
ries of high altitude balloon flights for direct detection of cosmic rays. The 
technique uses an initial charge-detecting layer for identifying the elemental 
charge, a target section for providing a nuclear interaction, and a calorimeter 
layer which collects the fraction of incident particle energy which appears as 
gamma-rays in the interaction (Si?^). In order to accumulate the required 
high energy statistics, JACEE has now flown emulsion chambers on 15 bal- 
loon flights (eight 1-2 day turnaround flights, two 5-6 day Australia-to-South 
America flights, and flve 9-15 day Antarctic circumpolar flights [7]). All but 
one of these have been successfully recovered. The total accumulated expo- 
sure is 1436 m^ hrs. The average flight altitude ranges from 3.5 to 5.5 g cm^. 
A single flight typically carries 2-6 emulsion blocks, each generally 40x50 
cm^. Fifty-eight emulsion blocks have been flown ; flfty-two have been recov- 
ered. The data discussed here come from JACEE flights 1-12, covering the 
results from the flrst 40 emulsion blocks (cumulative exposure 644 m^hrs). 
The total number of high energy events available for analysis from all flights 
is ~ 2 X 10"^. Of these, ~ 180 have energy exceeding 100 TeV per particle. 
The present analysis is based on 656 protons above 6 TeV and 414 helium 
nuclei above a total energy of 8 TeV per particle. The basic detector used in 
the JACEE experiments is a flnc grained emulsion [8] chamber typically con- 
taining approximately a hundred track sensitive nuclear emulsion plates and a 
three dimensional emulsion/x-ray/lead plate calorimeter. The lower part of the 
chamber is the calorimeter section, consisting of ~20 layers of emulsions and 
xray fllms interleaved with up to 8.5 radiation lengths of 1-2.5 mm thick lead 
plates. The calorimeter records singly charged particle tracks with a spatial 
resolution in the emulsion of better than 1 /im and individual photon cas- 
cades with a resolution of a few microns. High energy showers produce visible 
dark spots in the xray fllm, which are used to locate and trace the energetic 
cascades. On average, more than 400 events are detected per block with an 
optical density corresponding to a total energy in the electromagnetic shower 
HE^ > 1.5TeV for protons. In the original JACEE analyses [8-11] electron 
counts in the emulsion layers along the cascade were compared to a simulated 
shower development curve to determine the total electromagnetic energy de- 
posited in the calorimeter. In order to speed up the analysis for the large new 
data sample obtained from the long duration Antarctic flights, for these flights 
the energy is derived by fltting the observed darkness measured in the x-ray 
fllms to electromagnetic shower development curves to determine J^E^ [11,12]. 
The results of the measured fluxes of protons and helium nuclei into the knee 
region are shown in Figure 2. A flt to the unbinned integral spectra gives a 
proton spectrum with a differential power law index 2.80 ± 0.04 and a helium 
spectrum with a differential index 2.68 + 0.04/ — 0.06. There is no evidence of 
a break in either spectrum out to at least 40-90 TeV[14]. 
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Fig. 2. Direct measurements of cosmic ray components by the JACEE [13] and 
RUNJOB [15] experiments 



3.2 RUNJOB 



The RUssian-Nippon Joint Balloon collaboration (RUNJOB) has been car- 
rying out Trans-Siberian continental balloon flights since 1995, and eleven 
balloons have been launched from the Kamchatka peninsula, of which ten 
were successful. The payloads directly detect cosmic rays using a stack of 
passive detectors. Twenty emulsion chambers in total, each having an area 
of 40 cmx50 cm, have been recovered near the Volga region or west side of 
the Ural Mountains after an exposure of approximately 150 hours. The total 
exposure factor amounts to 575 m^ hrs. In this work the results based on 
the full analysis of the '95 and '96 experiments and one quarter of the '97 
experiment are given. The exposure factor used for the present data analysis 
amounts to 261.5 m^ hr. The performance of the balloon campaign and the 
data processing are presented elsewhere [16]. 

The payload mass of these flights is limited to 250 kg in each flight to ensure 
an altitude above ~10 g/cm^. The detector chambers consist of 4 parts: the 
Primary layer. Target layer. Spacer layer, and Calorimeter layer (see [15]). 
The primary layer is an emulsion used for charge identification, the calorime- 
ter consists of emulsion and x-ray films. The thickness of each layer is opti- 
mized based on the experience of our analysis. There is a wide spacer between 
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target module (made of lucite or stainless steel plate) and the calorimeter 
module, which allows a measurement of the opening angle of secondary 7-rays 
(tto 27) produced by the nuclear interaction at the target. The spacer thick- 
ness is typically 15-20 cm in the shower development direction, which gives a 
measurement of the opening angle for individual high energy 7-rays with en- 
ergy 20-50 TeV, taking the inclination effect into account. An estimate can be 
made of the total shower energy with ~100 TeV or more, where the elec- 
tromagnetic shower maximum may not be contained in the thin calorimeter, 
although this does depend on the inclination of the shower. The fraction of the 
incoming particle energy which is released into the 7-ray component can be 
estimated by using cascade shower theory, which is well-established through 
accelerator experiments [17]. The maximum darkness of the shower in the 
x-ray film is approximately proportional to the shower energy, Si?^, indepen- 
dent of the primary particle mass. The spot darkness is measured with use of 
a photo-densitometer, the details of which are summarized in [18] and [19]. 
However, the vertical thickness of the RUNJOB chamber is as thin as 4-5 
radiation lengths. A method to overcome this limitation has been discussed 
previously [20]. For heavy nuclei, the opening- angle method for fragments can 
be used. This method is summarized in detail in [21]. The results of the mea- 
sured fluxes of protons and helium into the knee region are shown in Figure 2. 



3.3 KASCADE 



The KASCADE (KArlsruhe Shower Core and Array DEtector) array is an 
air shower detector which consists of 252 detector stations in a 200 x 200m^ 
rectangular grid located near Karlsruhe, Germany [22]. The array comprises 
unshielded liquid scintillation detectors (e/7 detectors) and plastic scintillators 
as muon detectors below 10 cm of lead and 4 cm of steel. The total sensitive 
areas are 490 m^ for the e/7 detectors and 622 m^ for the muon detectors. 
In the center of the array a hadron calorimeter (16 x 20 m^ ) is built up, 
consisting of more than 40,000 liquid ionization chambers in 8 layers [23]. In 
between the calorimeter a trigger layer consisting of 456 scintillation detectors 
of 0.456 m^ each, measures the energy deposit and arrival times of muons with 
E>490 MeV. Below the calorimeter an array of position sensitive multiwire 
proportional chambers in two layers measures EAS muons with > 2.4GeV. 
Various methods and results of the analysis of KASCADE data arc described 
in detail in the proceedings of the last ICRC conference in Salt Lake City, 
Utah [24] or in [25], where the KASCADE contributions are sampled. The 
main observables of KASCADE per single shower are the so called shower 
sizes, i.e. total number of electrons and number of muons within the range 
of the core distance 40-200m, A''*'', local muon densities measured for the dif- 
ferent thresholds, and number and energy sum of reconstructed hadrons at 
the central detector. As a phenomenological result of KASCADE it should be 
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remarked that the frequency spectra of all these observables, i.e. for the dif- 
ferent particle components, show a clear kink at the same integral event rates. 
This is a strong hint for an astrophysical source of the knee phenomenon based 
on pure experimental data, since the same intensity of the flux corresponds 
to equal primary energy. The multi component measurements of KASCADE, 
especially the investigation of the hadronic component, allow the evaluation 
and improvement of hadronic interaction models in kinematical and energy 
ranges not covered by present collider experiments [26,27]. 

The shower sizes and N^J^ are used as input parameters for a neural net- 
work analysis to reconstruct the primary energy on an event by event basis. 
The necessary 'a priori' information in the form of probability density distri- 
butions are found by detailed Monte Carlo simulations with large statistics. 
The resulting energy spectrum depends on the high energy interaction model 
underlying the analyses. But the spectrum shows a clear kink at ~4-5 PeV 
and power law dependences below and above the knee. A simultaneous fit to 
the Ne and N^JI size spectra is performed for the reconstruction of the primary 
energy spectrum. The kernel function of this fit contains the size-energy cor- 
relations for two primary masses (proton and iron) obtained by Monte Carlo 
simulations. This approach leads to the all particle energy spectrum 
perposition of the spectra of light and heavy particles. Besides the use of global 
parameters like the shower sizes, sets of different parameters (describing dif- 
ferent shower particle components) are used for neural network and Bayesian 
decision analyses for showers with their axes within the central detector area. 
Examples of such observables are the number of reconstructed hadrons in the 
calorimeter {Eh >100 GeV), their reconstructed energy sum, the energy of the 
most energetic hadron ('leading particle' in the EAS), number of muons in the 
shower center {E^ >2 GeV), or parameters obtained by a fractal analysis of 
the hit pattern of muons and secondaries produced in the passive calorime- 
ter material. The latter ones are sensitive to the structure of the shower core 
which is mass sensitive due to different shower developments of light and heavy 
particles in the atmosphere. Monte Carlo statistics limit the number of pa- 
rameters which can be used for one multivariate analysis. Therefore a set of 
approaches using different observables are averaged in case of the actual result 
of the Bayesian analysis. The resulting classifications are corrected with mis- 
classification matrices leading to relative abundances. Afterwards the results 
are converted into distributions of the mean logarithmic mass. The resulting 
all-particle energy spectrum is shown in Figure 8 [29]. Figure 11 shows the 
mean logarithmic mass obtained by this approach versus primary energy [28] . 
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3.4 CASA-MIA 



The Chicago Air Shower Array - Michigan muon Array (CASA-MIA) de- 
tector [30] is an array of surface and underground plastic scintillators which 
measure the electromagnetic and muon components of air showers. The array 
is located on Dugway Proving Grounds, Utah, USA. The CASA-MIA data 
are used here to determine the probability, on an event by event basis, that a 
given data shower resulted from a 'light' or a 'heavy' primary. These distinc- 
tions are made by comparison to simulations of proton- or iron-induced air 
showers. A composition-independent measure of the energy is used to search 
for trends in the composition as a function of energy. 

In order to study composition, a detailed fit of all events is performed after a 
series of standard fits determines core distance, the electron and muon sizes 
and other shower quantities. Both the standard and detailed fits are based on 
the NKG and Greisen functions; the detailed fits allow more parameters to 
remain free. Three parameters sensitive to composition are extracted: the den- 
sity of surface particles Pe and the slope a of the lateral distribution near the 
core, and the density of muons at large core distance. The parameters pe, P/* 
and a are tabulated for data and for simulation events. Data events are classi- 
fied as 'iron-like' or 'proton-like' depending on whether these parameters most 
resemble iron or proton simulation events. The K Nearest Neighbor (KNN) 
test is employed to quantify this decision [31,32]. Each data event is placed in 
the three-parameter space defined by pe, and a. A large set of simulated 
iron and proton events also populate the space. Suitably normalized distances 
(computed in units of the variance of each parameter, with correlations in- 
cluded) between the datum and individual simulation events is calculated. 
The nuclear types of the five nearest simulation events are then tallied. The 
use of K=5 is optimal for this analysis [31]. Using too few simulation points is 
subject to fiuctuations, using too many would sample simulation points with 
very different primary energies. Each CASA-MIA event takes one of six pos- 
sible values, corresponding to whether it has 0-5 proton neighbors out of the 
five nearest neighbors examined. Using separate simulated calibration sets in 
place of data, it is found that more than 90% of events will have a majority of 
neighbors of their own species while about 50% have all nearest neighbors of 
their own kind. The 'proton resemblance' was normalized to the results of an- 
alyzing separate simulation calibration sets of pure iron and of pure protons. 
For example, near 100 TeV, the average fraction of proton nearest neighbors 
for a pure proton calibration sample is 68%; at the same energy a pure iron 
sample has only 16%. The average for CASA-MIA data at this energy is 51%. 
The computed proton resemblance is approximately proportional to (log(A)), 
the mean value of the logarithm of the atomic number of the sample [31,33]. 
These data are shown in Figure 11. 
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3.5 DICE 



The Dual Imaging Cherenkov Experiment (DICE) is a ground based air shower 
detector which is designed to have as little reliance as possible on the details of 
the air shower simulations and to have the capability of comparison with ex- 
isting direct measurements at 0.1 PeV to provide an assessment of the overall 
systematic error. Since the method of imaging hadronic showers in Cherenkov 
light is a relatively recent development, we provide some more detailed de- 
scription of the DICE detectors and operation. The two DICE telescopes are 
located at the CASA-MIA site in Dugway, Utah [30]). They each consist of a 
2m diameter f/1.16 spherical mirror with a focal plane detector of 256 close 
packed 40 mm hexagonal photomultiplicrs (PMTs) which provide ~ 1 deg pix- 
els in an overall field of view 16 deg x 13.5 deg centered about the vertical. The 
telescopes are on fixed mounts separated by 100 m (see [34,35]). 

Cosmic-ray events within the field of view produce a focal plane image at the 
PMT cluster which corresponds to the direction and intensity of Cherenkov 
light coming from the air shower. When the direction of the air shower core and 
the distance of the shower axis from the telescopes are known, simple geometry 
can be used to reconstruct the amount of light received from each altitude of 
the shower. The amount of Cherenkov light produced is strongly correlated 
with the number of electrons in the shower and is used to estimate the electron 
size as a function of depth in the atmosphere from which the location of shower 
maximum can be determined. This procedure is essentially geometrical and is 
independent of numerical simulations except for calculations which determine 
the angular distribution of Cherenkov light around the shower axis. For each 
air shower collected a simple time coincidence is used to identify the same 
event in both DICE clusters and CASA-MIA. Further requirements on the 
correlation of the DICE images with each other and with the CASA-MIA 
event geometry are used to reduce the overall probability of event mismatches 
between the detectors to ~10"^. The parameters for each shower are derived 
from these measured values. The accuracy of the shower core location derived 
by CASA is 1-3 m depending on the overall shower size. The measurement of 
the shower arrival direction is accurate to ~ 0.4 deg for larger showers with 
some degradation for lower energy events. 

Previous work with DICE estimated the shower energy by a simple translation 
from the total amount of Cherenkov light in the image and the geometry of the 
shower [36]. In this work [37] a more accurate estimate of energy is derived 
from a combination of the amount of Cherenkov light and the location of 
shower maximum size {X^^x-i in atmospheric depth g/cvo?) determined by each 
DICE telescope. This is desirable since the lateral distribution and intensity of 
Cherenkov light at a given total energy depends both on the primary particle 
mass, hence mean Xmaxi and the distance of the measurement from the shower 
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Fig. 3. Results from experiments reporting {Xmax) [39,37,53]. Also shown are the 
results from model calculations discussed in text. The SIBYLL [41] calculation is 
from Table 1, the QGSjet is from [42]. 



core. A fit for the total shower energy and primary particle mass is made to 
the geometry, Cherenkov size and X^ax location in the two DICE telescopes. 
The form of the Cherenkov size function used in these fits is derived from 
the results of simulations using the program CORSIKA 4.50 with the VENUS 
interaction model [38] . The derived X^ax fitting function has a constant shower 
elongation rate parameter and assumes simple superposition for providing 
primary particle mass dependence. 

With the acceptance cuts used during analysis the effective collection geometry 
is determined by the instrument Monte Carlo to be ~3300 m^ sr, making the 
overall collecting power ~125,000 m^ sr days. The results for the mean location 
of X^ax versus energy are shown in Figure 3. 
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3.6 BLANC A 



The Broad LAteral Non-imaging Cherenkov Array (BLANCA) takes advan- 
tage of the CASA-MIA particle array installation by augmenting it with 144 
angle-integrating Cherenkov detectors. Located in Dugway, Utah at an at- 
mospheric depth of 870 g cm~^, BLANCA uses the CASA trigger to collect 
Cherenkov light and records the Cherenkov lateral distribution from cosmic 
ray events in the energy range of the knee. The CASA trigger threshold im- 
poses an energy threshold of ~ 100 TeV on the Cherenkov array. However, 
BLANCA analysis uses events with a 200 TeV minimum to avoid composition 
bias introduced from the CASA trigger. 

Each BLANCA detector contains a large Winston cone [43] which concentrates 
the light striking an 880 cm^ entrance aperture onto a photomultiplier tube. 
The concentrator has a nominal half-angle of 12.5° and truncated length of 
60 cm. The Winston cones were aligned vertically with ~ 0.5° accuracy. A two- 
output preamplifier increases the dynamic range of the detector. The minimum 
detectable density of a typical BLANCA unit is approximately one blue photon 
per cm^. 

BLANCA operated between January 1997 and May 1998. After applying data 
quality cuts, about 460 hours of Cherenkov observations remain corresponding 
to an exposure of 212,000 m^ sr days. The detectors were calibrated on a 
nightly basis by using the circularly symmetric property of the cosmic ray air 
showers. This method was also verified with a portable LED fiasher. The LED 
system was then used in the laboratory to produce single photoelectrons in 
a PMT with good charge resolution to obtain the absolute calibration of the 
BLANCA system. A full description of the detector can be found in [44] . 

The lateral distribution of Cherenkov fight is fit to an empirical functional 
form to yield the normalization at 120 m (C120) and the exponential slope 

between 30 and 120 m (s). To interpret these measurements, a large set of air 
shower simulations were generated with the CORSIKA Monte Carlo package. 
To access the systematics inherent in the choice of interaction models, four 
separate hadronic models were used - QGSjet, VENUS, SIBYLL and HDPM. 
Ten thousand events were generated for each hadronic model in four primary 
groups - protons, helium, nitrogen and iron [45]. The simulated air shower 
events were then processed through a full BLANCA detector simulation [46] 
to obtain "fake" data. These fake data were used to find the optimum transfer 
functions between the observables C120 and s, and the energy and mass of the 
primary particle. The energy error distribution depends on the primary mass 
and energy. Assuming a mixed cosmic ray composition, the energy resolution 
for a single air shower is approximately 12% for a 200 TeV shower, falling to 
8% for energies above 5 PeV. In a similar way, a transfer function from slope, s. 
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to Xjnax was determined and the results are plotted in Figure 3. The quantity 
that is of direct astrophysical interest is the mean nuclear mass of the primary. 
In obtaining an absolute nuclear mass from the Cherenkov observables we 
chose to translate directly from (C120, s) to lii{A) to avoid compounding the 
errors that would occur by using primary energy and X^ax as intermediate 
quantities. For the purposes of this paper however, a translation through our 
Xmax values has been done for direct comparison with other groups. The 
results are shown in Figure 9. Full details of the BLANCA analysis are in [39]. 



3. 7 SPASE/VUL CAN/AMANDA 



An experiment to measure the electron, muon and air-Cherenkov components 
of ~1 PeV air-showers has been established at the geographic South Pole 
(688 g/cm^). The experiment comprises: the SPASE-2 scintillator array [47]; 
the VULCAN air-Cherenkov array [48] , and the deep under-ice muon detector 
AMANDA [49]. Simultaneous measurement of the electron size, high energy 
(>500 GeV) muon content and lateral distribution of Cherenkov light should 
allow a composition measurement that is relatively insensitive to model as- 
sumptions. 

SPASE-2 (a 1.6 X 10^ m^ array of plastic scintillator modules) provides an 
event trigger with a threshold of ~50 TcV for proton primaries. The SPASE- 
2 data are used to determine the shower core/direction to an rms accuracy 
of 4 m/l° at 1 PeV. Digitized waveforms from the nine element VULCAN 
array (similar in design to AIROBICC [50]) are used to reconstruct the lateral 
distribution of Cherenkov light. For a subset of events the shower axis passes 
through the AMANDA detector and the ice- Cherenkov signal can be used to 
estimate the > 500 GeV muon content of the shower. The aperture of the 
SPASE-2/ AMANDA-BIO telescope is falOO m^sr. 

The lateral distribution of Cherenkov light from air-showers is closely related 
to their longitudinal development [51] and can be used to provide an indi- 
rect measurement of the depth of shower maximum (Xmax) and primary en- 
ergy (E). The MOCCA air-shower simulation [52] (coupled with the SIBYLL 
hadronic interaction model [41]) has been used to deduce Xjnax{E) from the 
Cherenkov data obtained using the VULCAN array during 1997 and 1998 (see 
Figure 3). The systematic uncertainties in Xmax/E are ^ 15 g/cm^/20 %. This 
work is described in detail in [53] and [48]. A preliminary investigation of the 
distribution of X^ax at around 1 PeV suggests a ratio of (proton-|-helium) / all 
= O.Utom for the MOCCA/SIBYLL-17 model [54]. However the probabil- 
ity of compatibility with the data is low (2%). A much better fit is obtained 
by shifting the X^ax distributions to the mean values predicted by the COR- 
SIKA/VENUS model [55]: (proton+hehum)/all = 0.391^? with a probability 
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of 40%. 



For the subset of SPASE-2/VULCAN/AMANDA coincidences, a measure- 
ment of muon content from AMANDA can be combined with the VULCAN 
energy determination to estimate primary mass. Prehminary results from this 
method are consistent with the SPASE-2/ VULCAN X^ax result but have hm- 
ited statistics. These statistics will be increased by more than a factor of three 
by the addition of data from 1999 and 2000 in which the projected area of the 
AMANDA was increased. 

In addition the SPASE-2 and VULCAN surface arrays have been useful cal- 
ibration and surveying tools for the AMANDA detectors. It is anticipated 
that a surface air-shower detector would form a natural part of any future 
ice-Cherenkov neutrino observatory. 



4 Simulation Models for Air Showers 

The measurement of extensive air showers (EAS) is presently the only way to 
study cosmic rays with high statistics at energies near or above the knee of 
the all-particle spectrum. The properties of the primary cosmic rays have to 
be deduced from the shape and the particle content of the EAS. The energy 
of the primary is approximately reflected in the total number of secondaries 
produced and its mass can be deduced from the shower form. An iron primary 
e.g., consisting of 56 nucleons, would produce roughly a superposition of 56 
nucleon showers of energy E/56, thus being shorter with a shower maximum 
higher up in the atmosphere than a proton shower of the same energy. 

Since experiments at knee energies or above cannot be calibrated with a test 
beam the interpretation of EAS measurements is performed by comparison 
of the experimental data with model predictions of the shower development 
in the atmosphere. Therefore, unfortunately, quantitative results rely on the 
model assumptions and on the quality of the simulation of particle interactions 
and transport in the atmosphere. 

While the electromagnetic interaction (responsible for electromagnetic show- 
ers, ionization, Cherenkov light production, etc.) and the weak interaction 
(responsible for decays of unstable secondaries) are rather well understood, 
the major uncertainties in Monte Carlo simulations of EAS arise from the 
hadronic interaction models. The development of an EAS for a given en- 
ergy and primary is strongly dependent on the inelastic cross-sections Uinei 
of primary and secondary particles with air and on the fraction of the avail- 
able energy transferred into secondary particles (usually termed inelasticity). 
Large cross-sections and high inelasticity produce short showers, while small 
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cross-sections and low inelasticity produce long showers penetrating deep into 
the atmosphere. Unfortunately a systematic variation in cross-section and/or 
inelasticity in the model would directly translate into uncertainties of the pri- 
mary mass estimates obtained from the height of the shower maximum. 

The present theoretical understanding, in the form of QCD, cannot be applied 
to calculate the hadronic inelastic cross-section or the particle production from 
first principles. Therefore hadronic interaction models are usually a mixture 
of basic theoretical ideas and empirical parameterizations tuned to describe 
the experimental data at lower energies. Unfortunately none of the collider ex- 
periments register particles emitted into the extreme forward direction. These 
particles, being produced in collisions with rather small momentum transfer, 
carry most of the hadronic energy and arc of the greatest importance in the 
development of EAS, since they transport a large energy fraction deep into 
the atmosphere. Moreover, since cosmic ray energies exceed those accessible 
by accelerators by several orders of magnitude, models have to be extrapolated 
well beyond the range of firm knowledge. Here one has to rely on theoretical 
guidelines. 



4-1 CORSIKA and the hadronic interaction models 

The most widely used air shower model is the CORSIKA program [56,57]. 
It is a fully 4-dimensional Monte Carlo code, simulating in great detail the 
evolution of EAS in the atmosphere initiated by photons, nucleons, or nuclei. 
Single particles are tracked, accounting for energy loss, defiection due to mul- 
tiple scattering and the Earth's magnetic field, kinematically correct decays 
of unstable particles, and electromagnetic and hadronic interactions. COR- 
SIKA contains various hadronic interaction models which permits a compari- 
son and an estimate of systematic errors due to model assumptions and uncer- 
tainties in the hadronic interaction. The available models are HDPM [58,56], 
VENUS [59], QGSjct [60], SIBYLL [41], and DPMJET [61] and, for the ener- 
gies below Eiab = 80 GeV per nucleon, GHEISHA [62]. 

VENUS, QGSjet, and DPMJET model hadronic interactions of nucleons and 
nuclei on the basis of the Gribov-Regge theory (GRT) [63] , which most suc- 
cessfully describes elastic scattering and, via the optical theorem, the total 
hadronic cross-section as a function of energy. In GRT the observed rise of 
the cross-sections is a consequence of the exchange of multiple supercritical 
pomerons. Inelastic processes are described by cut pomerons leading to two 
color strings each, which fragment subsequently to color neutral hadrons. The 
probability of n exchanged and m cut pomerons is uniquely given by the the- 
ory. While this is common for all implementations of the GRT, they differ in 
the way the production and decay of the strings is realized. (A good introduc- 
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tion to GRT models is given in [59].) 



SIBYLL is a minijet model. It simulates a hadronic reaction as a combination 
of one underlying soft collision in which two strings are generated and a num- 
ber of minijet pairs leading to additional pairs of strings with higher p_[_ ends. 
In this model the rise of the cross-section is due to the minijet production 
only. 

HDPM is a purely phenomenological generator which uses detailed parameter- 
izations of pp collider data for particle production. The extension to reactions 
with nuclei and to energies beyond the collider energy range are somewhat 
arbitrary. 

GRT models simulate nucleus-nucleus coUisions in great detail allowing for 
multiple interactions of nucleons from projectile and target. The configurations 
of interacting projectile and target nucleons are determined with a Glauber- 
type geometrical picture of the nucleons in a nucleus [64]. HDPM and SIBYLL 
employ for nuclear projectiles the superposition model assuming an indepen- 
dent reaction for each of the projectile nucleon. 

QGSjct, DPMJET, and SIBYLL account for minijet production, which be- 
comes dominant at higher energies and leads to increasing cross-sections and 
more high-p^ particles. While in QGSjet and DPMJET minijets contribute 
only partially to the rise of the cross-section, they are solely responsible for it 
in SIBYLL. 

VENUS is the only model taking into account the interactions of intermediate 
strings and secondary hadrons with each other, leading to better agreement 
with final states as measured in collider experiments. In DPMJET these in- 
teractions are only implemented for secondary hadrons. 

For energies below 80 GeV GHEISHA is used, which was written around 
1985 to simulate the interactions of GeV secondaries from e+e~ coUisions 
with typical detector materials. GHEISHA is a phenomenological model, i.e. 
it parameterized particle production to agree with beam test results. 

Detailed comparisons of the models available in CORSIKA have been per- 
formed [55,65] reveahng differences on the 25% to 40% level and proving some 
of the models unable to describe aspects of the experimental results at all 
(e.g. SIBYLL and KASCADE results [66]). In the last three years some of the 
models were updated and new models were added to CORSIKA. 

QGSjet. 00 is an updated version of QGSjet with an increased Pomeron- nucleon 
coupling leading mainly to a cross-section with a 5% larger value at ^/s = 1800 
GeV, which is well within the variation of the experimental results [67-69] of 
about 10%. This was attempted to explain the trigger and hadron rates in 
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the KASCADE experiment [70]. However, while the EAS data could be repro- 
duced better by QGSjct.OO the rapidity distributions for p — p collisions were 
no longer reproduced by the model. Therefore it was decided to stick with the 
previous cross-sections. 

An error in the implementation of the DPMJET cross-sections was discovered 
and corrected which brought the values down by about 10%. 

A new version of SIBYLL (2.1) [71] was produced as a response to the deficit 
of muons in SIBYLL 1.6. SIBYLL now allows for multiple string production 
in soft collisions, which yields more muons and modifies the resulting cross- 
sections. Also better structure functions have been introduced to improve the 
extrapolations to high energies, and the treatment of nucleus-nucleus collisions 
was improved. Version 2.1 of SIBYLL is not yet officially released. The results 
shown here are based on a preliminary version which was available to the 
authors for evaluation. 

The model neXus [72] emerges as a new and more powerful program con- 
structed by the authors of VENUS and QGSjet. neXus is based on GRT 
but includes also hard interactions and attempts to describe all types of high 
energy interactions in a uniform way. neXus version 2 is implemented in the 
next CORSIKA release. 

URQMD [73] is an alternative to GHEISHA for low energy hadron interac- 
tions. It has been used in heavy ion collision experiments. It performs a more 
detailed modelling and would avoid problems with GHEISHA at a few GeV 
energies, but it turns out not to be free of problems either. 

4-2 Inelastic Cross-Sections 

A summary of experimental and predicted cross-sections is shown in Figure 4. 
Back in 1998 the variation of the p-air cross-section in experiments as well 
as in models amounted to about 25% at 3 x 10^^ eV, and to about 40% at 
10^*^ eV (shown in the left panel). For nucleus-air cross-sections the model 
uncertainty was only 10-15% due to the averaging effect of many nucleons 
inside the nucleus. The right panel shows the situation after some of the model 
authors have revised their cross-sections and new theoretical predictions have 
been published. HDPM, VENUS and MOCCA are not shown since they are 
no longer supported by their respective authors. 

Recent calculations of the p-air cross-sections from Frichtcr et al. [74] agree 
well with the QGSjet cross-sections. However, beyond 10^^ eV the predictions 
become slightly higher than those of QGSjet. Block et al. also [75] performed 
new calculations by using a QCD inspired parameterization of accelerator 
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Fig. 4. Inelastic p-air cross-sections, left: Situation in 1998 right: Situation in 2000 

data. They simultaneously fit the total p-p cross-sections, the ratio of real to 
imaginary part of the forward scattering amplitude p, and the nuclear slope 
parameter B. Then the p-p cross-sections were converted via Glauber theory 
into p-air cross-sections. Their result agrees very well with the cross-sections 
used in the QGSjet model. 

Extracting cross-sections directly from EAS measurements is not possible. 
Usually the attenuation length A of EAS in the atmosphere is measured which 
relates to the cross-sections crf^i^'^ via A = k - 14.6 mp/af^^" , where k depends 
strongly on the average inelasticity (and the inelasticity distribution) of p 
and TT reactions and on the energy. While early and naive models used the 
inelasticity as one of the free parameters, in CRT models the inelasticity dis- 
tribution emerges inevitably as a result of more fundamental properties of the 
interaction model. Usually k is assumed from an EAS model, which biases 
the final results for all the experimental cross-sections. In addition, different 
groups have been using different values for k. The experimental results can be 
easily brought into agreement with each other and with the theoretical calcu- 
lations by slightly modifying k [75]. For the AGASA and Fly's Eye data such 
a correction brought the experimental values in the energy range 10^^-10^^ eV 
down by 12-20%, which now are in good agreement with the cross-sections as 
predicted by Block et al. or QGSjet. 

As a result of various modifications the model predictions at knee energies 
agree now to within about 6%. 



4-3 Shower Development 



The cross-sections together with the inelasticity determine the longitudinal 
development of an EAS, which in turn is closely related to the most important 
shower observables: particle number at ground, height of shower maximum and 
total energy in the electromagnetic component. In Figure 5 average shower 
curves of 10^^ eV showers from different models are shown. 
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While for the old models the position of Xmax varied by about 7.5% for protons 
and 5.5% for iron showers, these number with the new models are only 1% 
and 0.5% respectively. Also the particle numbers at ground level agree now 
much better than two years ago. The variations for p and Fe showers were 
reduced from 80% and 50% to 14% and 3%, respectively. 

Calculations with URQMD instead of GHEISHA yield differences mainly in 
the low energy part of the showers. E.g. with URQMD there are about 30% 
more baryons in the energy range 2-10 GeV and a factor of 2 less for 0.1-2 
GeV than with GHEISHA. Also, while URQMD tends to produce 10% more 
muons below 4 GeV and 10% less above 4 GeV, the total number of muons 
stays approximately the same in both cases. 

The position of the shower maximum Xmax and its variation with energy are 
often used to estimate the mean mass and mass composition of primary cosmic 
rays. In Figure 6 calculated Xmax values are shown for QGSjet, SIBYLL 1.6 
and SIBYLL 2.1. 

None of the models show a linear relation as is often assumed in oversimplified 
arguments. As compared to SIBYLL 1.6, the new SIBYLL predicts about 20 
g/cm^ smaller Xmax values for both proton and iron induced showers. SIBYLL 
2.1 and QGSjet agree very well for iron showers and for proton showers until 
about 10^^ eV. Beyond this energy QGSjet gives smaller values. In the right 
panel of Figure 6 fluctuations are also shown. The shaded bands indicate 
the 1-a region around the average Xmax values. Note that the fluctuations of 
Xmax for proton showers are of the order of 100 to 60 g/cm^, i.e. not much 
smaller than the separation between protons and iron. These fluctuations make 
particle identification using Xmax for single events very difficult and require 
high statistics measurements to determine the average values reliably. In the 
case of QGSjet the differences between proton and iron showers decrease with 
energy, and consequently the bands for proton and iron start to overlap at 
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Fig. 5. Longitudinal shower development. The positions of the shower maxima are 
indicated at the upper rim of the picture, left: Situation in 1998 right: Situation in 
2000 
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Fig. 6. Height of shower nictximurn X^niax Si function of energy, left'. The curves 
for QGSjet, a prehminary version of SIBYLL 2.1 and SIBYLL 1.6 are shown, right: 
Average values together with the \-a band indicating the standard deviation of the 
Xmax distribution. [76] 

high energies, while for SIBYLL 2.1 they remain separate. Thus according to 
QGSjet the separation of different masses becomes more difficult with energy 
than according to SIBYLL 2.1. The X^ax values are summarised in Table 1. 



4-4 Collider Results 



Another, and sometimes underestimated, source of uncertainties are the errors 
of the collider results which are used to tune the interaction models. There 
are a few instances where the inherent precision of those results is revealed. 

The different inelastic p-p cross-sections as measured by experiments at Fer- 
milab were already mentioned [67-69]. Their results vary from 80.0 ± 2.2 mb 
to 71.7 ±2.0 mb at = 1800 GeV or Ejab = 1.7 x 10^^ eV. This corresponds 
to an 11% uncertainty already below the knee, and clearly much more for the 
highest energy cosmic rays. 

A second example is the pseudorapidity distribution of secondary particles 
produced in p-p collisions. The UA5 measurements from y/s = 200 to 900 
GeV [77] have been widely used to tune hadronic interaction models. GRT 
models always had severe difficulties reproducing the rapidity densities dN/d?7 
near 77 = and the slope in dN/dr] at around 77 = 4 at the same time. However, 
more recent measurements by Harr et al. [78] show a clearly flatter distribution 
than the UA5 ones, that fit very well to the GRT predictions (see Figure 7). 
In the near-forward region {rj ^ 4) the differences in the particle densities 
amount to about 25%. 
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Table 1 

Average X^ax and standard deviation for proton and iron showers and QGSjet and 
SIBYLL 2.1 models. The values are determined from 200 showers per energy point 
for iron and from 500 to 700 showers per energy for protons by fitting the shower 
maximum to each individual shower and averaging over the fit results [76]. 

4-5 Summary of air shower models 

In the last two years much has changed in EAS simulations. There is a clear 
trend of convergence between different hadronic interaction models, which, 
to some extent, is due to objective improvement of the physics put into the 
models. Presently QGSjet seems to be the best model showing a good overall 
agreement with a variety of experimental results. Typical differences in cross- 
sections or inelasiticities with other upcoming models (SIBYLL 2, neXus) are 
of the order of 10%. Unfortunately the number of independent models may 
be effectively reduced from 5 to 2 (neXus and SIBYLL 2), which will make 
it more difficult to estimate systematic uncertainties in the future, although a 
completely new version of DPMJET is in development [79] . Further improve- 
ment will be difficult since there is a variety of uncertainties at the 10% level in 
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Fig. 7. Pseudorapidity distributions as measured and simulated. 

collider results, in theoretical asumptions and, even larger ones, in experimen- 
tal air shower results, that enter into construction of EAS models. This makes 
it rather unlikely that EAS model predictions will get much better in the near 
future. A promising source of new information on hadronic interactions in the 
longer term will be the experiments at LHC and the new heavy ion colliders, 
even though CERN decided not to pursue the proposed experiments for very 
small angle secondaries. 

It is also unfortunate that many of the experiments investigating the compo- 
sition and spectra at the knee have ceased operation (CASA-MIA-BLANCA, 
DICE, EAS-TOP, HEGRA-AIROBICC). Also here it will be in future more 
difficult to combine different aspects of the problem and to estimate the sys- 
tematics. However, high statistics and good quality data are a good basis 
for substantial further improvement, perhaps with more model independent 
analysis techniques. 



5 Discussion of Results 

5.1 Energy Spectra 

Part of the intent of the workshop, and this report, is to examine how the 
various experiments are consistent in the measurement of composition in the 
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region of the knee. As a precursor to this comparison it is useful to see the 

measurements from the various experiments of the energy spectrum of cosmic 
rays in this region. These give some information on the systematic differences 
which may be present between the experiments, since the estabhshment of an 
absolute energy scale and the geometrical factor (m^sr s) of an experiment is 
often the most demanding task in the overall analysis effort. Figure 8 shows a 
comparison of the measured all-particle energy spectra from both direct and 
indirect measurements; only statistical errors are shown. The vertical intensity 
axis has been multiplied by E'^''^^, to emphasize differences from the overall 
spectrum of particles. The knee is clearly visible around 3 PeV in the air 
shower measurements. The overall flux normalization between the direct and 
indirect measurements is well within the errors, which are large and dominated 
by the limited geometrical factor of the direct measurements. Any of these 
measurements typically has a systematic error in the energy scale which is 
on the order 10-20%, this translates to a 20-40% shift in the overall intensity 
values normalization on this plot. The air shower measurements seem to agree 
almost too well with this level of expected systematic error. 



5.2 Composition estimates 

The derivation of composition from measurements such as those of mean 
shower maximum (Xmax), or derived from the electron and muon numbers, 
need a choice of the simulation to be used. As shown in Figure 3, both SIBYLL 
and QGSjet predict an overall logarithmic increase in the depth of shower 
maximum with energy. Here, to evaluate the size of systematic effects in the 
simulations, the SIBYLL data are taken from Table 1 [76] and the QGSjet 
data are from the recent work of Pryke [42] who has calculated similar quan- 
tities independently. The slope of this increase is fairly similar, but the overall 
normalization differences can have an impact on the mass estimates. 

As an example we can derive mean logarithmic mass composition from Fig- 
ure 3 by assuming the mass curves for other nuclei can be derived from a 
simple logarithmic interpolation between the proton and iron curves shown 
here. (This means that at any given energy {log{A)) = log{5Q) x {X^^^^ — 
{Xmax))/{XP^ax " ^max), where X^^^ and X^l^ are the model predictions 
and (Xjnax) is the experimentally measured value). Although this procedure is 
really only accurate if simple superposition is valid, it can give us an approx- 
imate basis on which to compare {log (A)) results between experiments. The 
results of this procedure for the (Xmax) results of Figure 3 using QGSjet are 
shown in Figure 9 where {loglO{A)) is shown versus energy, with the various 
data sets connected by dashed lines. Also shown here are the results from 
the direct measurements of the same quantity. There is good agreement be- 
tween the direct/indirect measurements, although the direct measurements 
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Fig. 8. All-particle energy spectra measurements near the 'knee' from various ex- 
periments presented at the workshop 

have poor statistics in the knee region. 

To evaluate the systematic errors in mass estimate, a similar procedure can be 
made with the SIBYLL model. The results of this are shown in Figure 10. Here 
the mean masses are larger, corresponding to the generally deeper shower max- 
ima expected for this model shown in Figure 3. However the general structure 
of the composition changes is preserved and again these are not inconsistent 
with the directly measured values. 

Within these Cherenkov air shower measurements, all compositions show a 
decrease in mean mass from around IPeV to 3PcV and two experiments 
(BLANCA and SPASE- VULCAN) show an increase in mean mass above 
3 PeV. A systematic shift in the energy scale of any of these air shower exper- 
iments will directly affect the inferred mass. The bar labelled 'systematic' at 
the left hand side of Figures 9 and 10 shows the shift in {log{A)) expected for a 
20% shift in energy scale. This is smaller than the systematic shift introduced 
between an anaylsis based on QGSjet or SIBYLL 2.1, the two curves shown 
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Fig. 9. Results on shower maximum analyzed with the QGSjet model 
in Figure 3. 

The goal of composition measurements near the knee is to look for apparent 
changes in the makeup of the all-particle spectrum which can reveal the origin 
of this feature. The most basic test which can be made is to examine whether 
the change in slope at the knee occurs at the same magnetic rigidity for all 
nuclear species. This would be expected if this feature has an origin either in 
an intrinsic steepening of the source spectrum of cosmic rays, possibly because 
diffusive shock acceleration becomes weaker, or because the particles begin to 
escape from the Galaxy more easily at these energies. An identification of the 
knee with a simple break in magnetic rigidity would significantly reduce the 
range of possible origins for this feature. 

The signature of a constant magnetic rigidity spectral bend is relatively easy 
to evaluate in terms of {loglO{A)). Since the all-particle spectrum of cos- 
mic rays measures particle intensities at similar total particle energies {E = 
'jAmc^, A=atomic mass, m=nucleon mass), at the same magnetic rigidity 
{R = ■jmcA/Z, since /3 ~ 1 at these energies) different nuclear species should 
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Fig. 10. Results on shower maximum analyzed with the SIBYLL model 

show a spectral bend at a total energy E — ZRc, where Z is the nuclear 

charge. For example, if a bend occurs for all nuclei at i?=lPV this will be ob- 
served at E=l PcV for protons and at E=2Q PeV for iron nuclei. The heavier 
nuclei are expected to bend at higher values of E, meaning that the mean 
mass should steadily increase across the knee region. 

An estimate of the size of this effect can be made by applying a simple model 
extrapolating from measured cosmic ray elemental abundances at lower ener- 
gies. An extension of a model of the Galaxy in which particles have an expo- 
nential distribution of pathlengths and an energy dependent leakage probabil- 
ity can be easily made into the 'knee' region [81]. These models produce an 
equilibrium intensity of particles in the Galaxy which are the result of subject- 
ing a source rigidity spectrum to losses by either spallation in the interstellar 
medium or escape from the Galactic confinement volume. Introducing a spec- 
tral break from oc i?"^-^ to oc R~^-^ at IPV produces a increase in (ZogrlO(A)) 
shown by the dot-dash curve labeled model in Figure 9 and Figure 10. 

This is a relatively robust prediction since the only assumption made is that 
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all nuclear species have the same source magnetic rigidity spectra. The source 
elemental abundances have been fit to lower energy data [81]. Any model 
which produces a bend in the rigidity spectra of this size will produce a similar 
increase in {loglO{A)) . The simple model does not agree well with either of 
the sets of derived mean compositions from the data. The measured values are 
consistently lighter than expected in two of the experiments. The Cherenkov 
experiments show a significant tendency to become lighter into the knee region 
at 3 PeV for all hadronic models. This does not agree with a simple rigidity 
bend origin for the knee. 

Figure 11 shows the estimates of {loglO{A)) derived from CASA-MIA using 
M0CCA/SIBYLL1.6 and KASCADE using CORSIKA/QGSjet. Also shown 
are the results from the direct measurements. The simple rigidity bend model is 
the dot-dashed curve discussed above. The derived compositions agree in gen- 
eral with the direct measurement normahzation. The KASCADE data show 
a trend towards greater mean mass even though there is a dip towards lighter 
values around 3 PeV similar to that seen with the Cherenkov experiments in 
Figure 9. CASA-MIA provides the most significant shift to heavy composition 
of all the experiments discussed here. 



5.3 Multi-Species Fits and Multi-Parameter Correlations 

The discussion in the last section is in terms of the mean value of a param- 
eter, be it {loglQ{A)) or {X^ax}- In fact, although fluctuations are severe, 
even a simple air shower experiment provides much more information in the 
form of the distributions and correlations of the observable quantities. This 
information can be used to improve the quality of the results — for instance 
providing a fractional composition breakdown into several elemental groups 
versus energy, rather than simply a mean mass. However, and probably more 
importantly given the highly indirect nature of the EAS experiments, the in- 
formation also provides a means to check in detail the consistency of the data 
with simulations, and hence to search for problems in the shower and detector 
Monte Carlos. 

In the past the hmited statistics of the data often meant that distribution 

analyses were not attempted — this is no longer the case. Additionally full 
shower and detector simulations, complete with all relevant fluctuations, were 
not carried out — this is now a standard part of any serious experiment. These 
two conditions being fuUfilled it becomes possible to compare the measured 
distribution of a composition sensitive parameter to the expected distributions 
for a range of primary masses. We normally assume that the range of possible 
primary particles in the knee region is bounded by protons at the light end 
and iron nuclei at the heavy. If this is the case the measured distribution must 
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Fig. 11. Mass measurements from KASCADE and CAS A- MIA compared with direct 
measurements. 

be reproducible by a sum of simulated distributions for proton, iron and some 
number of intermediate atomic masses. 

It is sometimes claimed that the severe inherent fluctuations of the air shower 
process limit the number of mass groups that may be fitted. In fact only the 
statistics of the data and the accuracy of the simulations impose such limits 
— in principle an arbitrary number of components can be deconvolved. In 
practice the impact of potential systematic errors in the simulations must be 
carefully assessed. 

If the observed distribution cannot be fit by an appropriate sum of simulated 
components then either the detector simulation is not adequately reproducing 
the measurement process, or the shower simulation is fiawed. If the shower 
simulation is suspected then, as discussed in Section 4, by far the most likely 
culprit is the high energy hadronic interaction model. In this way it can become 
possible to both choose an appropriate interaction model, and measure the 
primary composition simultaneously. 
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Fig. 12. Two examples of multi-species fits to air shower observables; BLANCA[39] 
at left, KASCADE[80] at right. See text for details. 

Two examples of sucessful multi species fits are shown in Figure 12. Both plots 
show the distribution of a shower observable which is stongly correlated with 
primary mass for a group of events in a band of primary energy around 3 PeV. 
At left we have the distribution of the exponential slope of the Cherenkov 
lateral distribution from BLANCA, and at right the ratio of the truncated 
muon number to electron number from KASCADE. The binned data is fit 
to a sum of the expected distributions for four primary species as simulated 
by CORSIKA-QGSjet, and after processing through the detector simulation, 
and event reconstruction procedures. The vertical axis is logarithmic in both 
cases to emphasize the tails of the distribution where inadequate simulations 
are most likely to be revealed. We can see that both experiments can be well 
fit by a four species model, and that neither appears to demand a sub-proton 
or trans-iron component. In fact both data sets demand a four species fit — 
the measured distribution cannot be adequately modeled with fewer (e.g. a 
simple proton/iron combination). It must be emphasized that such plots are 
a recent development, and a major step forward in air shower studies. 

A trend in recent experiments is to measure multiple parameters of the air 
shower on an event-by-event basis. It has sometimes been claimed that this can 
improve the resolution with which primary mass can be measured. However 
shower fluctuations are so severe that they easily dominate over instrumental 
resolution for any given parameter in a well designed experiment. Addition- 
ally the fluctuations are highly correlated between shower observables — in 
most cases what one is effectively measuring is the fluctuation of the elasticity 
(and multiplicity) of the first interaction. The real benefit of multi-parameter 
measurements is that not only the distributions of the indvidual observables, 
but also their cross correlations, can be used to compare with simulations. 
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Fig. 13. Two examples of event-by-event correlations between mass estimates near 
3 PeV derived from different air shower observables; BLANC A [39] at left from 
Cherenkov Lateral Distribution and Electron Muon sizes, DICE[37] at right from 
the Cherenkov image and Muon sizes. See text for details. 

Figure 13 shows two examples from the BLANCA and DICE extensions to 
the CASA-MIA system (see sections 3.6 and 3.5). BLANCA and DICE respec- 
tively measure the lateral and longitudinal distributions of Cherenkov light; 
CASA and MIA simultaneously measure the muon and electron sizes. Both 
plots are for a band of primary energy around 3 PeV (energy being estimated 
from the Cherenkov observables). The density of points in the xy plane is 
shown after making simple parametric transformations from observables to 
the logarithm of the primary mass [37,39]. Simulation plots are very similar 
for a mixed primary composition. The appearance of points on the figure with 
severly non-physical atomic masses serves to illustrate that the fluctuations 
of the air shower process are as large as the difference between the extreme 
ends of the mass scale which we are trying to measure. Note that the appar- 
ently better correlations in the case of BLANCA is at least partly due to the 
use of Nf./N^ rather than N^^ alone as the mass estimator — most of what is 
observed is that events with flat Cherenkov lateral distributions have larger 
electron sizes at ground level. 



6 Conclusions 

The experimental situation in measurements of cosmic-ray composition near 
the knee of the spectrum has beneflted in recent years from several new, rel- 
atively large, experiments. In the past these investigations were often un- 
dertaken with a few detectors and the simulation power for air showers was 
limited. Modern experiments have many more channels of information and 
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simulation codes have become extensive and sophisticated. In general the ex- 
periments discussed here do not suffer from statistical errors of sampling the 
shower, many measurements are taken for each shower. However the fact re- 
mains that trying to determine the identity of the nucleus which initiates an 
air shower is a very difficult task. A fundamental limit for this task seems 
to be connected to the accuracy of the air shower models which are used to 
derive the composition from the data. At present the systematic errors in es- 
timating the absolute energy scale and shifts between the air shower models 
limit the accuracy in {loglO{A)) to ~ 0.1. This is to be compared with an ex- 
pected shift of ~0.25 for a simple rigidity bend. There are mixed results from 
the experiments reported here for the existence of a simple bend. While the 
ground array measurements show a tendency for a heavier composition with 
increasing energy, the air cherenkov experiments show a tendency to lighten 
into the knee region with some increase in mean mass beyond 3 PeV. 

Several actions could improve this situation in the future. The accuracy of 
the modelling techniques can be improved if better accelerator results are 
available in the future as discussed previously. This also includes good mea- 
surements for small angle secondaries which at present do not exist. The direct 
measurements could be improved by longer exposures in space of high energy 
cosmic ray experiments. In particular this will help with the normalization of 
systematic errors in energy scale and mean mass measurements with air show- 
ers. Plans for using the Space Station for large payloads, such as the ACCESS 
mission which was discussed at the workshop (see http:/ /access.uchicago.edu), 
can improve the statistics near IPeV by about a factor of 50 over those pub- 
lished here. There are also other possible techniques to measure samples of 
cosmic rays at high energies. One of these discussed at the meeting was the 
prospect for measuring the direct emission of Cherenkov light produced by the 
particle before it interacts in the atmosphere using imaging techniques [82]. 
This could be a way to get significantly better elemental resolution at high 
energies with good statistics. Another promising approach is the multi compo- 
nent unfolding analysis discussed in Section 5.3. The K AS CADE experiment 
expects to extend this techniquJ^ to energies of 10^'^eV. 

The workshop would like to thank the Adler Planetarium in Chicago for all 
their help with organizing this meeting and providing an excellent venue for 
this event. Additional information discussed at the workshop but not directly 
presented here can be found at [http: / /knee. uchicago.edu, 



^ Since this workshop this work has been discussed at the 27th International Cosmic 
Ray Conference in Hamburg, Germany, 2001. Please see the conference proceedings 
for details [6]. 
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